Inverse swelling of a hydrophobic polymer in aqueous solution 
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We address the problem of inverse polymer swelling. This phenomenon, in which a collapsed 
polymer chain swells upon decreasing temperature, can be observed experimentally in so-called 
thermoreversible homopolymers in aqueous solution, and is believed to be related to the role of 
hydrophobicity in protein folding. We consider a lattice-fluid model of water, defined on a body- 
centered cubic lattice, which has been previously shown to account for most thermodynamic anoma- 
lies of water and of hydrophobic solvation for monomeric solutes. We represent the polymer as a 
self-avoiding walk on the same lattice, and investigate the resulting model at a first order approx- 
imation level, equivalent to the exact calculation on a Husimi lattice. Depending on interaction 
parameters and applied pressure, the model exhibits first and/or second order swelling transitions 
upon decreasing temperature. 

PACS numbers: 61.25.Hq, 05.50.+q, 87.10. +c 



I. INTRODUCTION 

From the experimental point of view, water is known 
to exhibit several thermodynamic anomalies, both as a 
pure substance 0, Q and as a solvent, in particular for 
non-polar (hydrophobic) chemical species 0, 0, Hj The 
transfer process of a small non-polar solute molecule in 
water is characterized by a positive solvation Gibbs free 
energy (it is thermodynamically unfavored), a negative 
solvation enthalpy (it is energetically favored), a nega- 
tive solvation entropy (it has an ordering effect), and 
a large positive solvation heat capacity [6j. More pre- 
cisely, for prototype hydrophobic species (that is, for in- 
stance, noble gases), solvation entropies and enthalpies, 
which are negative at room temperature, increase upon 
increasing temperature, and eventually become positive. 
Accordingly, the solvation Gibbs free energy displays a 
maximum as a function of temperature. These proper- 
ties define the so-called hydrophobic effect. 

It is quite well established that the hydrophobic effect 
is an important driving force for several biophysical pro- 
cesses , and in particular for protein folding j8| . This is 
the reason why it has attracted a high degree of attention 
in the last years, but a unified theoretical framework for 
this phenomenon does not exist yet. It has been observed 
that the native folded state of proteins is maximally sta- 
ble in the range of temperatures of living organisms, 
whereas it tends to be destabilized both by increasing and 
by decreasing temperature. For globular proteins, it is 
eventually possible to observe complete denaturation also 
upon decreasing temperature [9|, [lfj . This phenomenon 
is denoted as cold unfolding, and has a simple analogue 
in so-called thermoreversible homopolymers, which ex- 
hibit a transition from a collapsed to a swollen state, 
upon decreasing temperature 0, 0, 0. We denote 
this kind of phase transition as inverse swelling. All this 
phenomenology is qualitatively consistent with the pre- 
viously mentioned existence of a maximum in the sta- 
bilizing force, i.e., hydrophobic repulsion, as a function 



of temperature. Such reasoning is of course unrigorous, 
since the hydrophobic effect is not a real interaction, but 
an effective repulsion, resulting from an average over the 
degrees of freedom of water. 

Cold unfolding and inverse swelling have attracted 
some attention from the theoretical point of view 0, 
III 13 IH 13 E El El El IH, and they have been 
a starting point for highlighting the importance of wa- 
ter structure details for modeling protein folding |24| . 
One possible way of investigation consists of computer 
simulations, based on very detailed (all-atom) interac- 
tion potentials. Unfortunately, simulations are generally 
limited by the large computational effort required. In- 
vestigation of a full protein model with explicit water 
is still out of reach [25j , and also for homopolymers the 
analysis is generally limited to relatively short chains or 
small parameter regions [2(| 1271 128| . A complementary 
approach involves investigations of simplified models, ei- 
ther on- or off-lattice. Several models of this kind have 
been proposed ElElEElEEIlSlllllllll. 
As previously mentioned, it turns out that, in order to 
reproduce inverse swelling (or cold unfolding), it is nec- 
essary to take into account, even in extremely simplified 
fashions, the degrees of freedom of water. Nevertheless, 
such descriptions are generally proposed ad-hoc for this 
kind of problem. The focus is mainly on the polymer, 
whereas water models by themselves would not be satis- 
factory models of water. 

Conversely, in this article, we start by considering a 
lattice-fluid model of water, whose thermodynamic prop- 
erties have been previously investigated in detail |23.l30| . 
This model predicts most of the thermodynamic anoma- 
lies of real water at constant pressure (a temperature of 
maximum density, a minimum of isothermal compress- 
ibility and specific heat), and also a liquid- liquid phase 
separation in the supercooled region, and a second criti- 
cal point |29j . Moreover, a corresponding model for aque- 
ous solutions of ideally inert monomeric solutes turns out 
to exhibit the above mentioned fingerprints of hydropho- 
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bicity, and in particular the maximum in the solvation 
free energy as a function of temperature [3(j. We de- 
scribe the polymer in solution as a lattice self-avoiding 
walk, whose steps connect nearest neighbor sites. We 
assume that each visited site, representing a monomcric 
unit, interacts with water and with other non-consecutive 
monomeric units, in the same way as the elementary so- 
lute of the the original model does. In some sense, our 
model is unbiased. 

The water model is a simplified version of the one 
proposed by Roberts and Debenedetti 0, 113, HI]. It 
is defined on the body-centered cubic lattice, and wa- 
ter molecules possess four equivalent bonding arms ar- 
ranged in a tetrahedral symmetry. According to this 
model, the microscopic description of water anomalies is 
based on the competition between an isotropic (van der 
Waals-like) interaction and a highly directional (hydro- 
gen bonding) interaction, and on the difference between 
the respective optimal interaction distances. In the lat- 
tice environment, the latter is taken into account by a 
trick, that is, a weakening of hydrogen bonds by neigh- 
boring water molecules. 

The previously mentioned results on this model have 
been obtained by a first order approximation on a tetra- 
hedral cluster (e quiv alent to the exact calculation for 
a Husimi lattice |3J| made up of tetrahedral building 
blocks), but turn out to compare quite well with Monte 
Carlo simulations (2^- The same kind of approxima- 
tion has been independently verified to be rather accu- 
rate for the case of a semiflexible self-interacting polymer 
chain [3^ with no explicit solvent. We thus expect that 
the calculation presented here, based on the same ap- 
proximation technique, could produce reliable results as 
well, although it is being applied to a more complicated 
model involving water-polymer interactions. 



II. THE MODEL 

Let us now describe the model in detail. As previously 
mentioned, it is defined on the bec lattice. Each site 
may be empty or occupied by a water molecule (w) or 
by a polymer segment (s). Since we find it more conve- 
nient to work in the grand-canonical ensemble, each of 
the two species contributes to the hamiltonian by a dif- 
ferent chemical potential term — /i x , where x S {w, s}. An 
attractive (Van der Waals) energy — e xy < is assigned 
to any pair of nearest neighbor (NN) sites occupied by 
x, y molecules, where x and y run over {w,s}. Of course, 
the interaction energy between polymer segments — e ss 
is taken into account only if the segments are not con- 
secutive along the chain. Water molecules possess four 
equivalent arms that can form hydrogen bonds, arranged 
in a tetrahedral symmetry, so that they can point to- 
wards 4 out of 8 NNs of a given site. A hydrogen bond 
is formed whenever two NN molecules have a bonding 
arm pointing to each other, yielding an energy —rj < 0. 
Only 2 different water configurations can form hydrogen 




FIG. 1: Husimi lattice structure and bonding water config- 
urations i = 1 (left) and i = 2 (right). Views from top are 
meant to explain configuration symbols used in Table [i] 

bonds (see Fig. QJ; w more configurations are allowed, 
in which water molecules cannot form bonds (the w pa- 
rameter is related to the bond-breaking entropy). An 
energy increase ijc x /6, with c x £ [0, 1], is added for each 
of the 6 sites closest to a bond (i.e., 3 out of 6 next 
nearest neighbors of each bonded molecule) , occupied by 
an x molecule. As far as water molecules are concerned, 
the weakening parameter mainly accounts for the fact 
that hydrogen bonds are most favorably formed when 
water molecules are located at a certain distance, larger 
than the optimal Van der Waals distance. Therefore, 
if too many molecules are present, the average distance 
among them decreases, and hydrogen bonds are (on av- 
erage) weakened. Moreover, the presence of an external 
molecule may perturb the electronic density, resulting in 
a lowered hydrogen bond strength as well. A weaken- 
ing parameter c s for polymer segments could take into 
account possible perturbation effects for a generic chem- 
ical species. Nevertheless, in this paper we shall consider 
only the case of an ideal monomer with c s = 0. 

The bec lattice model can be conveniently replaced by 
an analogous model defined on a Husimi lattice made 
up of irregular tetrahedral building blocks (see Fig. QJ, 
as explained in detail in Refs. l29l and I30L The latter 
model, which can be studied exactly, is equivalent to a 
generalized first order approximation (on the tetrahedral 
cluster) for the original model |34| . From now on, we shall 
refer to the Husimi lattice model only. The hamiltonian 
of the system can be written as a sum over the tetrahedral 
clusters 

W= ^ Ti ioi i i^i s , (1) 

{ot,0,l,S) 

where Tiijki is the contribution of each cluster, which 
will be referred to as tetrahedron hamiltonian, and the 
subscripts i a , ip, i 7 , ig label site configurations for the 4 
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vertices (a, f3,j,S). It is understood that the latter are 
enumerated in a conventional order, for instance clock- 
wise, with reference to Fig. Q (bottom). Accordingly, 
(a,/3), (/?, 7), (7,<S), and (S, a) turn out to be NN pairs, 
whereas (a, 7) and (/3, S) are next nearest neighbors. 

Site configurations are reported in Table [IJ where 
i = 1,2,3 denote water molecule configurations (defined 
as in Refs. l29tl3(J) and £ = 4, — ,8 denote local configura- 
tions of the polymer chain. As far as polymer configura- 
tions are concerned, let us notice that they are more con- 
veniently defined with respect to the tetrahedron, than 
with respect to a fixed reference frame. This trick allows 
us to write a unique form for the tetrahedron hamilto- 
nian, which otherwise would be orientation-dependent. 
Moreover, although in the isotropy hypothesis every local 
chain configuration should be equally probable, we need 
to distinguish different classes of configurations (each one 
with a given multiplicity) in order to impose connectiv- 
ity constraints, as it will be clarified below. As a con- 
sequence, the Husimi lattice approximation introduces a 
small artifact, in that the probability of the configura- 
tions i = 4, 6 (in which the chain locally stays in the 
same cluster) turns out to be slightly different from the 
probability of the configurations i = 5, 7, 8 (in which the 
chain passes from one cluster to another). We expect 
that such artifact is negligible, according to the results 
of Ref.lM 

The tetrahedron hamiltonian can be written as 
7~iijki — Hijki + Hjku + Hkuj + Hujk + Lijki, (2) 
where 



Hiju = -e xy n x ^n yj bfb- - rjhij 1 1 - c x x ' - — — J 
--fi x n x>i + l(bf®bj). (3) 

Let us analyze the various terms appearing in these equa- 
tions. 

The first line of Eq. J2J is basically equivalent to the 
corresponding term for the mixture model of water and 
a simple monomeric solute [see Eq. (3) in Ref. |3(J , and 
includes the Van der Waals and the hydrogen bonding 
terms. Occupation variables are defined as n Xj ,; = 1 
if the configuration i corresponds to a chemical species x, 
and n Xj i = otherwise (see Table , whereas hij are 
hydrogen bond variables, defined as hij — 1 if the pair 
configuration represents a hydrogen bond (i.e., if 

i = 1 and j = 2), and hij = otherwise. It is understood 
that repeated x and y indices are summed over their pos- 
sible values w, s. The only difference with respect to the 
monomeric solute case is related to the fact that, as pre- 
viously mentioned, we have to exclude Van der Waals 
interactions between consecutive polymer segments. To 
do so, we have defined the bond numbers bf (resp. b~[), 
which are boolean variables equal to 1 if the configura- 
tion i represents a polymer segment forming a chemical 
bond in the clockwise (resp. counterclockwise) direction 



TABLE I: Site configurations, with corresponding labels (i), 
multiplicities (u)i), occupation variables (n„ t i, ns,i), and chem- 
ical bond numbers (bf,b~). Polymer configurations are de- 
fined with reference to the cluster denoted by a dot; 
denote the configurations "viewed" by the other 3 clusters: 
i = empty site; i = 1,2 bonding water; i — 3 non-bonding 
water; i = 4 segments in the same cluster, out of the ref- 
erence one; i = 5 segments in different clusters, out of the 
reference one; i = 6 segments in the same cluster, the ref- 
erence one; i = 7 segments in different clusters, one in the 
reference cluster ("clockwise" direction); i — 8 segments in 
different clusters, one in the reference cluster ( "counterclock- 
wise" direction); 
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of the reference cluster, and otherwise (see Table QJ. 
In this way, the multiplying term bfbj (where the over- 
line denotes boolean negation) is if represents two 
consecutive segments, and 1 otherwise. 

The second line of Eq. © includes the chemical po- 
tential contributions (multiplied by 1/4 to avoid over- 
counting) and an infinite energy penalty / — > 00, as- 
signed to configurations which violate connectivity 
constraint, i.e., when i wants to form a chemical bond 
(in the clockwise direction) and j does not (in the coun- 
terclockwise direction), or vice versa. This is obtained 
by the exclusive-or factor (bf ®bj). Let us notice that 
the infinite energy penalty can be treated numerically, 
since we have to deal only with Boltzmann factors of the 
tetrahedron hamiltonian. Finally, the additive term Lijki 
in Eq. (J2J is another constraint term, needed to forbid 
short loops on tetrahedra. It is therefore simply defined 
as £6666 — 00 and L^ki — otherwise. 



III. THE FIRST ORDER (HUSIMI LATTICE) 
APPROXIMATION 

As mentioned in the Introduction, we perform the ex- 
act analysis of the Husimi lattice model, which is equiva- 
lent to the generalized first order approximation on tetra- 
hedral clusters. The calculation closely follows the one of 
the ordinary mixture model, performed in Ref.l3Ct so that 
we do not give much detail. Basically, one has to take 
into account that the system is locally tree-like, which 
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allows to write a recursion equation [3J], for so-called 
partial partition functions, denned as follows. Let us as- 
sume that the lattice is actually a tree (Husimi tree), 
and let us consider just one branch, with the correspond- 
ing partial hamiltonian, obtained by Eq. Q with the 
sum restricted to tetrahedra in the branch. The partial 
partition function Qi is a sum of Boltzmann weights of 
the partial hamiltonian, taken over the configurations of 
the branch except the base site (this is why the partial 
partition function depends on the base site configuration 
variable i). It is convenient to define normalized partial 
partition functions qi oc Qi, for instance in such a way 
that 



53* 

i=0 



(4) 



If the branch becomes infinite, i.e., in the thermody- 
namic limit, and in the homogeneity hypothesis, the 
subbranches attached to the first tetrahedral cluster are 
equivalent to the main one, so that one can write the 
recursion equation 



ii = y 1 E • 

J=0 



8 8 

^ Wk ^ wie 



fc=0 



Ti-ijki/T 



1=0 



Y[q^Qksqis, (5) 



where the sum runs over configuration variables in the 
tetrahedron except i, T is the temperature expressed in 
energy units, and y is a normalization constant. 

Let us notice a subtle but important difference with 
respect to the corresponding equation for the monomeric 
solute [Eq. (8) in Ref . l30l| . As mentioned in the previous 
section, we have defined local chain configurations with 
reference to a given cluster, which allows us to write one 
single tetrahedron hamiltonian, independently of orien- 
tation. Therefore, when we consider the operation of 
"attaching branches" to a base cluster, we have to take 
into account configurations of sites in the cluster as they 
are "viewed" by the attached branches. To do so, we 
have to write a product over £ = /, //, ///, which denote 
precisely the different views, according to Table [I] 

The recursion equation can be iterated numerically to 
determine a fixed point. The site configuration probabil- 
ities pi can be computed by considering the operation of 
attaching 4 branches to the given site, yielding 



Pi = Z 1 qiY\q i i, 
where the normalization constant is determined as 



(6) 



8 

£ 

i=0 



Wiqi \\ qa ■ 



(7) 



From the knowledge of site configuration probabilities pi, 
one can easily compute the densities, i.e., the probabili- 
ties p w , p s that a site is occupied by a water molecule or 



by a polymer segment, respectively. For x G {w,s}, we 
have 



^ Wiping, 

i=0 



(8) 



where the occupation numbers explicitly given in 

Tabled 

In the presence of multiple fixed points (which can be 
reached from different initial conditions), i.e., in the pres- 
ence of coexistence phenomena, the first order transition 
can be determined by minimizing the (grand-canonical) 
free energy per site 



u> = —T (In y — 2 In z) . 



(9) 



where y and z are the normalization constants of the re- 
cursion equation and of the site probability distribution, 
respectively. The derivation of this expression requires 
some manipulations, and is left to the Appendix. From 
the knowledge of the free energy, one can in principle 
compute all other thermodynamic properties. Assuming 
the volume per site equal to 1, pressure can be expressed 
in energy units as P — — u>. 



IV. RESULTS 

First of all, we have to fix a set of model parameters. 
As far as water is concerned, we choose the values of 
our previous investigation [2!j: e W w/^ = 0.25, w = 20, 
and c w = 0.5. The hydrogen bond energy r\ is taken as 
the energy unit. The water- water van der Waals energy 
e ww is significantly smaller. The multiplicity w of non- 
bonding water configurations is large, to mimic the high 
directionality of hydrogen bonds. 

With this set of parameters, the model provides a qual- 
itatively consistent description of the phase diagram and 
thermodynamic anomalies of pure water as already 
mentioned in the Introduction. For the liquid phase, one 
observes a density maximum as a function of tempera- 
ture at fixed pressure, i.e., a change of sign in the isobaric 
thermal expansion coefficient 



ap 



din p„ 



dT 



(10) 



The temperature of maximum density slightly decreases 
upon increasing pressure, as observed in experiments. 
The isothermal compressibility 



Kt 



dhiycv 



dP 



(11) 



and the isobaric specific heat exhibit a minimum as a 
function of temperature, at constant pressure. More- 
over, in the supercooled regime, the model predicts the 
so-called second critical point, which has been conjec- 
tured and observed in simulations 01 , and of which also 
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FIG. 2: Pressure (P/rj) vs temperature (T/rf) phase diagram 
for pure water. The solid line denotes the liquid- vapor transi- 
tion line, the thin solid line denotes the liquid phase spinodal, 
whereas the thin dotted line denotes the temperature of max- 
imum density locus. Parameter values: e ww /?7 = 0.25, w — 20, 
c w = 0.5. 



some experimental evidences have been found [361]. The 
liquid-vapor spinodal displays no reentrance in the posi- 
tive pressure half-plane, whereas the temperature of max- 
imum density locus exhibits a peculiar "nose-shaped" 
reentrance. Let us recall that this scenario qualitatively 
agrees with recent molecular dynamics simulations of 
water phenomenological potentials in the negative pres- 
sure region 2], whereas the reentrant spinodal scenario 
is an older conjecture invoked to explain thermodynamic 
anomalies of liquid water |37l l38j . The plausibility of the 
latter has been recently a matter of debate [3!J |4(| , but 
seems to be definitely ruled out. 

We summarize the results obtained by our model for 
pure water in Fig. |2| (phase diagram) and Fig. (response 
functions), in order to delimit the temperature range in 
which we ought to expect the (stable or metastable) liq- 
uid water regime. Let us notice that the response func- 
tions correctly exhibit a divergent behavior, upon ap- 
proaching the liquid phase spinodal. 

Upon inserting an ideal inert solute, i.e., adding a 
monomeric s chemical species, characterized by e ws = 
e ss = and c 5 = 0, we have shown that the model is able 
to reproduce in a qualitatively accurate way also the sol- 
vation properties of simple hydrophobic solutes, such as 
noble gases [3(| . Since in Ref. |3(] we have not considered 
the present set of water parameters, we report here the 
results, in terms of solvation free energy, entropy, and 
enthalphy as a function of temperature (Fig.^J. 

The solvation Gibbs free energy per molecule Ag* , i.e., 
the transfer free energy of a solute molecule from vapor 
to liquid phase at coexistence, can be computed as 



Ag* = -Tln-^y, 

Ps 



(12) 




0.28 0.30 0.32 0.34 0.36 0.38 0.40 0.42 

temperature 

FIG. 3: Pure water response functions at fixed pressure 
(P/rj — 0.05) as a function of temperature (T/rf). Upper 
graph: isobaric thermal expansion coefficient (r/ctp). Lower 
graph: isothermal compressibility (t]kt)- Parameter values: 
Eww/f) = 0.25, w = 20, c w = 0.5. 
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0.28 0.30 0.32 0.34 0.36 0.38 0.40 0.42 

temperature 

FIG. 4: Transfer energy functions (E/rf) vs temperature 
(T/rf) for a single monomer at liquid- vapor coexistence: E = 
Ag" (solid line), E = TAs* (dashed line), E = Ah* (dotted 
line). Parameter values: e„„/ri — 0.25, w = 20, c w = 0.5, 
£ws = e ss = 0, c s = 0. 



phases, respectively. The solvation entropy is defined as 

8Ag 



As* 



dT 



(13) 



where and ps denote solute densities in the two 



and the solvation enthalpy can then be computed as 

Ah* = Ag* +TAs*. (14) 

The typical experimental temperature range for several 
hydrophobic compounds is around 0°C to 300° C Ac- 
cording to pure water results, in our model this roughly 
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corresponds to the range T/n w 0.33 (just below the tem- 
perature of maximum density) to T/rj ~ 0.38 (about half 
way between the previous temperature and the critical 
temperature). In this range, the solvation Gibbs free en- 
ergy is positive and displays a broad maximum, whereas 
entropy and enthalpy are negative at low temperatures 
and positive at high temperatures, in agreement with ex- 
periments [6(. Upon approaching the critical point, the 
free energy tends to zero, whereas the entropy and en- 
thalphy diverge. 

For the polymeric solute, we carry out the investigation 
as discussed in the previous Section. We compute, as a 
composition variable, the segment molar fraction, defined 
as 



(15) 



We perform the analysis at constant pressure P/r] = 0.05, 
i.e., the same value for which we have computed pure wa- 
ter response functions. At fixed temperature, we always 
observe a transition between a phase with x ~ (pure 
water) at lower values of the segment chemical potential 
/i s , and another phase with a certain positive fraction 
x of polymer segments {polymerized phase) at higher /i 5 
values. In order to determine this transition, we have 
programmed a numerical procedure, which adjusts the 
water chemical potential in order to fix the pressure 
of the two phases, and then finds the /z s transition value, 
at which the water chemical potentials of the two phases 
are equal. The polymerized phase at the transition de- 
scribes the behavior of an isolated polymer chain in solu- 
tion . Accordingly, the x value at this point is an indi- 
cator of the polymer chain compactness, as a function of 
temperature. If x = 0, the transition is continuous, and 
the polymer is in a completely swollen state, whereas, if 
x > 0, the transition is first order, and the polymer is in 
a more or less collapsed state. 

As far as interaction parameters are concerned, we 
have first investigated a completely inert (hydrophobic) 
polymer, assuming e ws = e ss = and c 5 = 0, and then we 
have turned on a slight attractive water-segment interac- 
tion (e W5 > 0). Depending on the value of this parameter, 
we have obtained different interesting behaviors. The re- 
sults are summarized in Figs. |S] and |H| 

From e ws = to e W5 /n w 0.04, we observe a first- 
order inverse swelling, since the segment molar fraction x 
abruptly jumps from a high value down to zero (Fig. [SJ- 
In the ]i s vs T phase diagram (Fig. top panel), this 
transition corresponds to a critical end-point, which lies 
among two different polymerized phases with different 
densities, and the water phase. Upon decreasing the seg- 
ment chemical potential /x s , the low temperature poly- 
merized phase undergoes a second order transition to the 
pure water phase, corresponding to the swollen state. On 
the contrary, the high temperature polymerized phase un- 
dergoes a first order transition to the pure water phase, 
corresponding to the collapsed state. For e W5 — 0, the 
inverse swelling transition occurs around T/rj « 0.227. 




0.28 0.30 0.32 0.34 0.36 0.38 0.40 0.42 

temperature 

FIG. 5: Polymer compactness (segment molar fraction x for 
the polymerized phase at the transition with the pure water 
phase) as a function of temperature (T/n) at fixed pressure 
(P/r] = 0.05), for different values of e ws /rj (reported on each 
plot). Other parameter values: e ww /r] = 0.25, w — 20, c w = 
0.5, e ss = 0, c s = 0. 
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FIG. 6: Segment chemical potential (fis/n) vs temperature 
(T/n) phase diagrams at fixed pressure (P/r] = 0.05), for 
different values of e ws /r] (reported on each graph). Other pa- 
rameter values: e„„/?7 = 0.25, w = 20, c„ = 0.5, e ss = 0, 
c s = 0. Solid and dashed lines denote first and second order 
transitions, respectively. 



At such a low temperature, we expect that the liquid 
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water phase predicted by our model is actually unstable, 
so that the observed phase behavior is to be viewed as 
an extrapolation j^. Upon increasing 6 W s, the inverse 
swelling transition temperature increases. It is interest- 
ing to notice that also an "ordinary" swelling transition 
can be observed at high temperature. Such transition 
coincides with the liquid-vapor transition of pure water, 
and it is therefore completely unaffected by the water- 
segment interaction parameter e ws . In our opinion, such 
feature is a clear evidence that the polymer collapse is in- 
duced by the presence of a dense (liquid) aqueous solvent 
around it, i.e., by the hydrophobic effect. 

For e W5 /?7 > 0.04, the behavior changes, in that in- 
verse swelling occurs via an intermediate, moderately 
collapsed phase (Fig. 0. Upon decreasing temperature, 
we first observe a discontinuous jump of the segment 
molar fraction x from a high value to a smaller value, 
and subsequently a continuous transition to a completely 
swollen state with x = 0. In the fi s vs T phase diagram 
(Fig. El middle panel), the former transition corresponds 
to a triple point among the two polymerized phases and 
the water phase. The latter transition corresponds to a 
"tricritical" point, at which the transition between the 
low temperature polymerized phases and the pure water 
phase changes from first to second order. Such behavior 
resembles a 8 point |4l| in which the role of tempera- 
ture is reversed. Upon increasing e ws , the temperatures 
of both transitions increase, and the amplitude of the 
discontinuous one is progressively reduced. 

Around e ws /?7 « 0.10, we observe a totally continu- 
ous inverse swelling, but a reminiscence of the first order 
jump can still be observed (Fig. EJ. Such a behavior cor- 
responds, in the fj, s vs T phase diagram (Fig. [SI bottom 
panel), to the fact that the transition between the two 
different polymerized phases ends in a critical point be- 
fore encountering the transition line with the pure water 
phase. For e ws /r) w 0.11 and above, we can no longer 
observe a collapsed phase, and the polymer is completely 
swollen at all temperatures. Accordingly, in the fj, s vs T 
phase diagram, the transition between the polymerized 
phase and the pure water phase is second order at all 
temperatures. 

It is noticeable that the effect of pressure on the poly- 
mer behavior is somehow similar to the effect of the 
water-segment attractive energy e ws , as shown in Fig.0 
The first-order swelling transition temperature increases 
upon increasing pressure. The transition becomes less 
and less abrupt, and is finally smoothed out at high pres- 
sure values, i.e., in the pressure vs temperature phase di- 
agram, the transition line ends in a critical point. Let us 
notice, by the way, that a similar critical behavior has 
been predicted by a simplified model of cold unfolding 
for proteins [lfij . 

We can also observe some differences. For instance, 
the high temperature swelling transition, which was unaf- 
fected by the water-segment interaction at fixed pressure, 
increases upon increasing pressure, although, at pressure 
values lower than the critical one, it still coincides with 



0.8-, 




CD 



CO 

0.4- 
0.2 -. 



0.28 0.30 0.32 0.34 0.36 0.38 0.40 0.42 

temperature 

FIG. 7: Polymer compactness (segment molar fraction x for 
the polymerized phase at the transition with the pure water 
phase) as a function of temperature {T/rj), for different pres- 
sure values P/rj (reported on each plot). Parameter values: 
£ww/?7 = 0.25, w = 20, c w = 0.5, e„ s /?7 = 0.07, e ss = 0, c s = 0. 



the liquid-vapor transition. This means that indeed, at 
high temperature, pressure tends to make the polymer 
coil more compact. The same holds true for the low tem- 
perature, moderately collapsed phase. This effect is re- 
versed at very high pressures. 

It is important to stress that the different polymer be- 
haviors described above partially occur in a temperature 
region which is likely to be unreachable by experiments. 
According to Fig. [21 the response functions of pure water 
display peaks, but in actual experiments only their high 
temperature "sides" can be observed. The real peaks are 
believed to lie below the homogeneous nucleation temper- 
ature By the way, this is the reason why there has 
been a debate about the possibility of a reentrant spin- 
odal for liquid water (in the latter case the peaks would 
become actual divergences). As far as the polymer prob- 
lem is concerned, we conclude that only temperatures 
above T/rj ss 0.32 can be experimentally investigated. 
We shall return to this issue in the following. 

As previously mentioned, the present model provides 
a qualitatively good description of hydrophobic solvation 
thermodynamics for a simple monomeric solute. In par- 
ticular, if the monomer is ideally non- interacting, the 
transfer free energy displays a maximum as a function 
of temperature (see Fig. 0J. This is expected to be a 
key ingredient for inverse swelling of a hydrophobic poly- 
mer. Nevertheless, in the previous investigation, we have 
taken into account the effect of a weak water-monomer 
interaction. Such assumption allows us to observe inverse 
swelling at experimentally accessible temperatures. For 
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FIG. 8: Transfer free energy (Ag* /rj) as a function of temper- 
ature (T/rj) for a single monomer at liquid-vapor coexistence, 
for e ws /r? = 0, 0.03, 0.05, 0.08, 0.09, 0.10. Other parameter val- 
ues: e-mj/rj — 0.25, w = 20, c„ = 0.5, e ss = 0, c s = 0. 



completeness, we then compute the solvation free ener- 
gies for the monomeric solute, for all the employed values 
of e ws /ry. The results are reported in Fig. [HJ 

As a general trend, the curves are shifted towards lower 
energy values, upon increasing e W5 , while the maximum 
becomes broader. Such behavior qualitatively accounts 
for progressive hindering of the collapsed state observed 
in the polymer model. Nevertheless, we find it remark- 
able that, from the solvation free energy alone, one would 
not at all expect the previously described variety of be- 
haviors. In our opinion, this is a striking evidence of the 
fact that inverse swelling cannot be satisfactorily mod- 
eled by effective interaction potentials, meant to take 
into account a temperature-dependent "quality of the sol- 
vent" , but it is rather to be ascribed to a subtler interplay 
between polymer and water degrees of freedom. 



V. CONCLUSIONS 

In this paper, we have proposed and investigated a lat- 
tice model of a linear hydrophobic homopolymer in aque- 
ous solution. As far as water is concerned, the model 
is a simplified version of a model proposed by Roberts 
and Debenedetti [3l| , in which we have removed distinc- 
tion between donors and acceptors, which is generally be- 
lieved not to play a crucial role in the thermodynamics of 
hydrogen bonding With respect to other simplified 
model addressed to study the same kind of systems and 
phenomenology, the main new ingredient of the present 
work is that we have extended a previously studied model 
of water and aqueous solutions to the case in which the 
solute is a polymer, without other ad-hoc assumptions. 

We observe the inverse swelling phenomenon, i.e., a 
transition from a collapsed to a swollen state upon de- 
creasing temperature. For an ideally inert (hydrophobic) 
polymer, such transition is first order, and occurs at very 
low temperatures, where we expect that liquid water is 



unstable. In this respect, our model results can be viewed 
just as an extrapolation. Nevertheless, the transition can 
be moved into the temperature region of liquid water by 
an attractive water-monomer interaction. In these condi- 
tions, inverse swelling is no longer purely first order, but 
occurs via a discontinuous transition to a moderately col- 
lapsed phase, followed by a O-like transition to a com- 
pletely swollen state. The latter transition likely stays 
in the unobservable region, so that we expect that only 
first order swelling could be observed. Upon increasing 
pressure, the first order transition temperature increases, 
the transition is smoothed, and eventually disappears at 
high pressure values. Let us briefly discuss the relevance 
of such findings to experimental results. 

First of all, the typical system in which the inverse 
swelling phenomenon can be observed is an aqueous so- 
lution of a thermoreversible homopoly mer, such as poly- 
N-isopropylacrylamide (PNIPAM) [Hill 03- lt is in- 
teresting to notice that, for PNIPAM, inverse swelling oc- 
curs actually via a first order transition at a temperature 
close to 27° C, at ordinary pressures. In the framework 
of our model, a similar behavior can be reproduced, as 
previously mentioned, only in the presence of a water- 
segment attractive interaction. For PNIPAM, the exis- 
tence of such attraction can be justified by the fact that 
the monomeric unit (NIPAM) contains in the side chain a 
polar (hydrophilic) peptide group, beside ten non-polar 
(hydrophobic) hydrogens (the structure is reported for 
instance in Ref. [11). It is noteworthy that the assump- 
tion of such kind of interaction had to be made as well 
in a previous phcnomcnological treatment of the inverse 
swelling phenomenon for PNIPAM [2^. In the cited ar- 
ticle, the solvation free energy for NIPAM was assumed 
as an effective temperature-dependent interaction with 
water, and, as a consequence, inverse swelling turned out 
to be a 0-like point, i.e., a continuous transition. Con- 
versely, from our investigation it clearly turns out that 
the information contained in the solvation free energy is 
insufficient to discriminate the type of transition. 

A purely hydrophobic short polymer chain in water 
has been recently investigated in great detail, by means 
of molecular dynamics simulations |2^. Beside several 
structural results, it has turned out that inverse swelling 
cannot be observed at ordinary pressures, but only at 
very high pressures. The transition is abrupt enough to 
suppose that it becomes discontinuous in the thermo- 
dynamic limit, and the transition temperature increases 
with pressure. These results qualitatively agree with our 
findings, although for our model the effect of pressure 
is not sufficient to move the transition in the observable 
region, for the completely inert polymer. 

We have mentioned in the Introduction that some in- 
terest on the inverse swelling phenomenon has been at- 
tracted by the fact that globular proteins undergo a 
somehow similar phenomenon, known as cold unfolding. 
Our model is definitely too simple to describe the physics 
of proteins, first of all since a protein is a heteropolymer, 
in which the monomeric units have different degrees of 
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affinity with water. Nonetheless, some interesting analo- 
gies exist. In most cases, the cold unfolding transition 
of proteins can be observed in the temperature regime 
of stable liquid water, in the presence of denaturants 
(urea) [9|, or b y ap plying high pressures (pressure de- 
naturation) E5- Extrapolations to zero denaturant 
concentration, or to ambient pressure, suggest a transi- 
tion temperature below the freezing temperature of wa- 
ter. Both effects seem to be somehow consistent with our 
results. On the one hand, a denaturant could be roughly 
viewed as something which on average makes the protein 
more hydrophilic, i.e., a water-segment attractive interac- 
tion. On the other hand, we have observed that pressure 
tends to increase the transition temperature, although 
this effect is relatively weak. 

Let us finally recall that the results presented here have 
been obtained by an approximate technique, i.e., a gen- 
eralized first-order approximation. Such kind of calcula- 
tion requires quite a negligible computational effort, but 
usually provides qualitatively reliable results. In partic- 
ular, the accuracy of the approximation has been previ- 
ously verified, by a comparison with Monte Carlo sim- 
ulations, both for the pure water model [2^| and, in- 
dependently, for an ordinary self-avoiding walk model 
with stiffness [^J . Nevertheless, a similar test for the 
full model would be welcome, but this goes beyond the 
scope of the present paper, and is left to future work. 

APPENDIX: EQUILIBRIUM FREE ENERGY 

Hereafter, we justify the formula employed in the text 
to evaluate the equilibrium (grand-canonical) free energy 
density. It is known that the exact solution of a Husimi 
lattice can be obtained also in a variational approach, 
as a minimum of a generalized Bethe free energy |8 it . 
Eq. @ could be written in principle from the general 
free energy formula reported in Ref . l34l but this formula 
has been proved only for the case of scalar configura- 
tion variables. Here, although with some tricks we can 
label configurations with scalar quantities, polymer con- 
figurations actually behaves as vector variables, which 
is revealed by the presence of "primed" configurations 
in Eqs. 10 and 0. Therefore, we have to extend the 
proof, restarting from the variational Bethe free energy. 
For our particular system, in the homogeneity hypothe- 
sis, the Bethe free energy per site can be written as 

uj/T = (Hijki/T + lnP ijkl ) - 3 (ln Pi ) , (A.l) 

where (■) denotes an ensemble average over configuration 
variables, and Piju is the tetrahedron probability distri- 
bution. Let us notice that the previous formula holds for 
a generic free energy density, not necessarily the grand- 



canonical one. The tetrahedron probability distribution 
can be written as a function of the partial partition func- 
tions in the following form 

/// 

P m = 2- 1 e- w «"/ r JJffit^-e^g,*, (A.2) 
f=' 

which results by considering the operation of attaching 
3 branches to each site of the given tetrahedron. The 
normalization constant is 

8 8 8 8 /// 

Z = ^ Wi ^2 w j E Wk E w i e ~ Hl3kl /T II IV- ?j« 1k( <7;« • 

i=0 j=0 fc=0 i=o £=/ 

(A.3) 

Replacing Eqs. (|A.2(1 and © into HA.1|) . and taking into 
account the linearity of the average operation, one ob- 
tains by simple algebra 

uj/T= -lnZ + 31nz + $, (A.4) 

where 

in I > in 8 

* = E ( ln — ) = E E w &* ( ln «* - ln ?*) • ( A - 5 ) 

£=/ \ * / £=/ i=0 

Expanding the sum, and inserting the multiplicity values 
Wi reported in Table [IJ one can write 

$ = 3 (p e - pi) ln g 4 + 3 (p 4 - p 6 ) In q 6 

+12(p 7 +p 8 -2p 5 )lng 5 (A.6) 
+6 (2p 5 +p 8 - 3p 7 )lnq 7 + 6 (2p 5 + p 7 - 3p$) \nq 8 , 

where the terms associated to configurations i = 0, 1, 2, 3 
(empty site or water) disappear, because in this case i = 
i' = i" = i'" holds (see Table [IJ . As mentioned in the 
text, one can observe that p^ = pe and ps = p 7 = ps, 
whence $ = 0. Moreover, making use of Eqs. (|A.3|) . J5j, 
and J7J), it is possible to show that 

Z = yz, (A.7) 

where y is the normalization constant of the recursion 
equation (J5J, which has to be computed at each itera- 
tion, and z is given by Eq. Q. Replacing the previous 
equation into Eq. (|A.4|I with $ = 0, one finally obtains 
Eq. 10, which we aimed to prove. 
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